General setup

Setup chunk

Load libraries

knitr::opts_chunk$set(fig.width = 8)
knitr::opts_knit$set(root.dir = normalizePath(".."))
knitr::opts_knit$get("root.dir")
[1] "/nas/groups/treutlein/USERS/tomasgomes/projects/pallium_evo"

Set colours for cell types and regions

library(Seurat)
library(harmony)
library(ggplot2)
library(Matrix)
library(mgcv)
library(foreach)
library(doParallel)
library(parallel)
library(URD)
library(slingshot)
library(RColorBrewer)

Load data

meta = read.csv("data/annotations/axolotl_all_umeta.csv", 
                header = T, row.names = 1)
cols_cc = c(
#epen
"#12400c", "#2d6624","#1d4f15", "#174711", "#2d6624", "#3d7f33", "#3b7b30", "#468b3b", "#4f9843","#5dae50", "#66bb58", "#72cd64", "#306a26", "#78d669", "#81e472",
#gaba
"#700209", "#75090e","#7a0f13", "#801517", "#851a1b", "#8a1f1f", "#902423", "#952927", "#9a2d2c","#a03230", "#a53634", "#aa3a39", "#b03f3d","#b54342", "#ba4846", "#c04c4b", "#c5504f", "#ca5554", "#d05959", "#d55e5e","#73050c", "#780c11","#8d2221", "#982b2a","#a23432", "#a83837", "#b2413f", "#b84544", "#bd4a49", "#c85352", #"#cd5756",
#glut
"#054674", "#134d7b","#1d5481", "#265a88", "#2e618e", "#73a4cb", "#366995", "#3e709c", "#4677a2","#4d7ea9", "#5586b0", "#5c8db7", "#6495bd","#6b9cc4", "#7bacd2", "#8ebfe4", "#96c7eb", "#9ecff2", "#18507e", "#18507e","#2a5e8b", "#497ba6","#5889b3", "#6fa0c8","#7fafd6", "#6091ba", "#5182ac", "#3a6c98", "#a6d7f9",
#npc
"#ffb120", "#feb72a","#fdbc34", "#fcc13d", "#fbc745", "#facc4e", "#f9d156", "#f8d65f", "#f8da68","#f7df70", "#f7e479", "#f7e882", "#f7ed8a", "#f7f193", "#eca319"
)
ccnames = unique(sort(meta$cellclusters))
names(cols_cc) = c(ccnames[grepl("epen", ccnames)], ccnames[grepl("GABA", ccnames)],ccnames[grepl("glut", ccnames)],ccnames[grepl("npc", ccnames)])

reg_cols = c("other/unknown_pred" = "#C7CCC7", 
             "medial" = "#52168D", "medial_pred" = "#661CB0", 
             "dorsal" = "#C56007", "dorsal_pred" = "#ED7307", 
             "lateral" = "#118392", "lateral_pred" = "#16A3B6")
reg_cols_simp = c("medial" = "#52168D", "dorsal" = "#C56007", "lateral" = "#118392")

colors <- colorRampPalette(brewer.pal(11,'Spectral')[-6])(100)

Calculate pseudotime

Process data

Define groups

div_srat = readRDS("data/expression/axolotl_reclust/Edu_1_2_4_6_8_12_fil_highvarfeat.RDS")
pred_meta = read.csv("results/Div-seq/divseq_predicted_metadata.csv", 
                     header = T, row.names = 1)[,4:5]
colnames(pred_meta) = c("region", "cellclusters")
div_srat = AddMetaData(div_srat, metadata = pred_meta)

Subset data

common_cells = c("epen_clus_3", "epen_clus_4")
r_cells = list("hippocampus" = c("glut_SUBSET_0", "glut_SUBSET_4", "glut_SUBSET_5",
                                 "glut_SUBSET_7", "npc_SUBSET_7","npc_SUBSET_11",
                                 "glut_SUBSET_14", "glut_SUBSET_15", "glut_SUBSET_12",
                                 "glut_SUBSET_16", "glut_SUBSET_17", "glut_SUBSET_3",
                                 "npc_SUBSET_1", "npc_SUBSET_3", "npc_SUBSET_9"),
    "lc" = c("glut_SUBSET_9","npc_SUBSET_7","npc_SUBSET_4", "glut_SUBSET_21", "glut_SUBSET_2",
             "glut_SUBSET_25","npc_SUBSET_14","npc_SUBSET_0"),
    "cl111" = c("glut_SUBSET_1", "glut_SUBSET_11", "npc_SUBSET_2","npc_SUBSET_4",
                "npc_SUBSET_7", "npc_SUBSET_13"),
    "eomes" = c("npc_SUBSET_7","npc_SUBSET_4", "glut_SUBSET_10","glut_SUBSET_22"),
    "cl8620_ep" = c("glut_SUBSET_8","glut_SUBSET_6","glut_SUBSET_20", "epen_clus_1",
                    "epen_clus_7", "epen_clus_14"))

Quick processing

sub_srat = list()
for(n in names(r_cells)){
  sub_srat[[n]] = div_srat[,div_srat$cellclusters %in% c(common_cells, r_cells[[n]])]
}

URD

Prepare URD objects

sub_urd = list()
for(n in names(sub_srat)){
  sub_urd[[n]] = sub_srat[[n]]
  sub_urd[[n]] = createURD(count.data = sub_urd[[n]]@assays$RNA@counts, 
                           meta = sub_urd[[n]]@meta.data, min.cells=3, min.counts=3)
  sub_urd[[n]]@group.ids$cellclusters = sub_urd[[n]]@meta$cellclusters
  sub_urd[[n]]@group.ids$subclasses = sub_urd[[n]]@meta$subclasses
}

Calculate DRs

for(n in names(sub_urd)){
  sub_urd[[n]]@var.genes = findVariableGenes(sub_urd[[n]], set.object.var.genes=F,do.plot=F,
                                           diffCV.cutoff=0.3, mean.min=.005, mean.max=100)

  sub_urd[[n]] = calcPCA(sub_urd[[n]], mp.factor = 2)
  pcSDPlot(sub_urd[[n]])
  sub_urd[[n]] = calcTsne(object = sub_urd[[n]])
  plotDim(sub_urd[[n]], "cellclusters")
  
  sub_urd[[n]] = calcDM(sub_urd[[n]], knn = 150, sigma=16)
  
  plotDimArray(sub_urd[[n]], reduction.use = "dm", dims.to.plot = 1:8, 
             outer.title = "Diffusion Map (Sigma 16, 150 NNs): Stage", label="cellclusters",
             plot.title="", legend=F)

  plotDim(sub_urd[[n]], "cellclusters", transitions.plot = 10000)
}

Pseudotime calculation

for(n in names(sub_urd)){
  sub_urd[[n]]@var.genes = findVariableGenes(sub_urd[[n]], set.object.var.genes=F,do.plot=F,
                                           diffCV.cutoff=0.3, mean.min=.005, mean.max=100)

  sub_urd[[n]] = calcPCA(sub_urd[[n]], mp.factor = 2)
  pcSDPlot(sub_urd[[n]])
  sub_urd[[n]] = calcTsne(object = sub_urd[[n]])
  plotDim(sub_urd[[n]], "cellclusters")
  
  sub_urd[[n]] = calcDM(sub_urd[[n]], knn = 150, sigma=16)
  
  plotDimArray(sub_urd[[n]], reduction.use = "dm", dims.to.plot = 1:8, 
             outer.title = "Diffusion Map (Sigma 16, 150 NNs): Stage", label="cellclusters",
             plot.title="", legend=F)

  plotDim(sub_urd[[n]], "cellclusters", transitions.plot = 10000)
}
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
[1] "2022-05-25 16:34:24: Centering and scaling data."
[1] "2022-05-25 16:34:33: Removing genes with no variation."
[1] "2022-05-25 16:34:35: Calculating PCA."
[1] "2022-05-25 16:49:52: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 3.09454943131581"
[1] "14 PCs have eigenvalues larger than 2 times null upper bound."
[1] "Storing 28 PCs."
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
[1] "Using provided global sigma 16"
Warning in DiffusionMap(data.use, sigma = sigma.use, k = knn, n_eigs = dcs.store,  :
  You have 4087 genes. Consider passing e.g. n_pcs = 50 to speed up computation.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
[1] "2022-05-25 17:06:56: Centering and scaling data."
[1] "2022-05-25 17:07:01: Removing genes with no variation."
[1] "2022-05-25 17:07:02: Calculating PCA."
[1] "2022-05-25 17:17:34: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 3.46600896608338"
[1] "14 PCs have eigenvalues larger than 2 times null upper bound."
[1] "Storing 28 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
[1] "Using provided global sigma 16"
Warning in DiffusionMap(data.use, sigma = sigma.use, k = knn, n_eigs = dcs.store,  :
  You have 3920 genes. Consider passing e.g. n_pcs = 50 to speed up computation.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
[1] "2022-05-25 17:26:42: Centering and scaling data."
[1] "2022-05-25 17:26:47: Removing genes with no variation."
[1] "2022-05-25 17:26:49: Calculating PCA."
[1] "2022-05-25 17:36:30: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 3.33792950730749"
[1] "13 PCs have eigenvalues larger than 2 times null upper bound."
[1] "Storing 26 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
[1] "Using provided global sigma 16"
Warning in DiffusionMap(data.use, sigma = sigma.use, k = knn, n_eigs = dcs.store,  :
  You have 3711 genes. Consider passing e.g. n_pcs = 50 to speed up computation.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
[1] "2022-05-25 17:45:46: Centering and scaling data."
[1] "2022-05-25 17:45:51: Removing genes with no variation."
[1] "2022-05-25 17:45:52: Calculating PCA."
[1] "2022-05-25 17:54:20: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 3.68127706582674"
[1] "12 PCs have eigenvalues larger than 2 times null upper bound."
[1] "Storing 24 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
[1] "Using provided global sigma 16"
Warning in DiffusionMap(data.use, sigma = sigma.use, k = knn, n_eigs = dcs.store,  :
  You have 3764 genes. Consider passing e.g. n_pcs = 50 to speed up computation.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning in densfun(x, parm[1], parm[2], ...) : NaNs produced
[1] "2022-05-25 18:00:38: Centering and scaling data."
[1] "2022-05-25 18:00:42: Removing genes with no variation."
[1] "2022-05-25 18:00:45: Calculating PCA."
[1] "2022-05-25 18:06:45: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 3.92503020569977"
[1] "12 PCs have eigenvalues larger than 2 times null upper bound."
[1] "Storing 24 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
[1] "Using provided global sigma 16"
Warning in DiffusionMap(data.use, sigma = sigma.use, k = knn, n_eigs = dcs.store,  :
  You have 3690 genes. Consider passing e.g. n_pcs = 50 to speed up computation.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

Save

end_cc = list("hippocampus" = c("glut_SUBSET_0", "glut_SUBSET_7"),
              "lc" = c("glut_SUBSET_9","glut_SUBSET_2"),
              "cl111" = c("glut_SUBSET_1", "glut_SUBSET_11"),
              "eomes" = c("glut_SUBSET_10","glut_SUBSET_22"),
              "cl8620_ep" = c("glut_SUBSET_8","glut_SUBSET_6"))
subset_dat_l = list()
floods_4_l = list()
for(n in names(sub_urd)){
  print(n)
  # Here we use all cells from the first stage as the root
  root.cells = cellsInCluster(sub_urd[[n]], clustering = "cellclusters", "epen_clus_4")
  
  # Then we run 'flood' simulations
  floods = floodPseudotime(sub_urd[[n]], root.cells = root.cells, n=250,
                           minimum.cells.flooded = 2, verbose=F)
  # fix excessive NA
  floods = floods[,colSums(!is.na(floods))>1000]
  
  # The we process the simulations into a pseudotime
  sub_urd[[n]] = floodPseudotimeProcess(sub_urd[[n]], floods, floods.name="pseudotime_4")
  
  pseudotimePlotStabilityOverall(sub_urd[[n]])
  
  plotDists(sub_urd[[n]], "pseudotime_4", "cellclusters", plot.title="Pseudotime by stage")
  
  # Create a subsetted object of just those cells from the final stage
  subsetdat = urdSubset(sub_urd[[n]], 
                        cells.keep=cellsInCluster(sub_urd[[n]], "cellclusters", end_cc[[n]]))
  subset_dat_l[[n]] = subsetdat
  floods_4_l[[n]] = floods
}
[1] "hippocampus"
[1] "lc"
[1] "cl111"
[1] "eomes"
[1] "cl8620_ep"

Determine new identities for end tips

save(subset_dat_l, sub_urd, floods_4_l, file = "data/processed/URD/URD_Div_lists.RData")

Random walk simulations

axials_4_l = list()
for(n in names(sub_urd)){
  # Use the variable genes that were calculated only on the final group of stages (which
  # contain the last stage).
  subset_dat_l[[n]]@var.genes = sub_urd[[n]]@var.genes[sub_urd[[n]]@var.genes %in% rownames(subset_dat_l[[n]]@logupx.data)]
  
  # Calculate PCA and tSNE
  subset_dat_l[[n]] = calcPCA(subset_dat_l[[n]], mp.factor = 1.5)
  pcSDPlot(subset_dat_l[[n]])
  
  set.seed(20)
  subsetdat = calcTsne(subset_dat_l[[n]])
  
  # Calculate graph clustering of these cells
  subset_dat_l[[n]] = graphClustering(subset_dat_l[[n]], 
                                      num.nn = 50, do.jaccard=T, method="Louvain")
  plotDim(subset_dat_l[[n]], "Louvain-50", 
          plot.title = "Louvain (50 NN) graph clustering", point.size=3)
  plotDim(subset_dat_l[[n]], "cellclusters", plot.title = "cellclusters", point.size=3)
  
  
  # Copy cluster identities from axial.6somite object to a new clustering ("tip.clusters") in the full axial object.
  sub_urd[[n]]@group.ids[rownames(subset_dat_l[[n]]@group.ids), "tip.clusters"] = subset_dat_l[[n]]@group.ids$`Louvain-50`
  
  # Determine the parameters of the logistic used to bias the transition probabilities. 
  # The procedure is relatively robust to this parameter, but the cell numbers may need to be
  # modified for larger or smaller data sets.
  axial.ptlogistic = pseudotimeDetermineLogistic(sub_urd[[n]], "pseudotime_4",
                                                 optimal.cells.forward=20, max.cells.back=20,
                                                 do.plot = T)
  
  axial.biased.tm = as.matrix(pseudotimeWeightTransitionMatrix(sub_urd[[n]], "pseudotime_4",
                                                               logistic.params=axial.ptlogistic))
  
  axials_4_l[[n]] = list("log" = axial.ptlogistic, "tm" = axial.biased.tm)
}
[1] "2022-05-26 00:16:23: Centering and scaling data."
[1] "2022-05-26 00:16:24: Removing genes with no variation."
[1] "2022-05-26 00:16:24: Calculating PCA."
[1] "2022-05-26 00:16:38: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 9.51302270671831"
[1] "4 PCs have eigenvalues larger than 1.5 times null upper bound."
[1] "Storing 8 PCs."
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

[1] "Mean pseudotime back (~20 cells) 0.00161580965872476"
[1] "Chance of accepted move to equal pseudotime is 0.5"
[1] "Mean pseudotime forward (~20 cells) -0.00161580965872476"
[1] "2022-05-26 00:16:51: Centering and scaling data."
[1] "2022-05-26 00:16:52: Removing genes with no variation."
[1] "2022-05-26 00:16:52: Calculating PCA."
[1] "2022-05-26 00:17:03: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 9.67601122396447"
[1] "5 PCs have eigenvalues larger than 1.5 times null upper bound."
[1] "Storing 10 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

[1] "Mean pseudotime back (~20 cells) 0.00215910967243115"
[1] "Chance of accepted move to equal pseudotime is 0.5"
[1] "Mean pseudotime forward (~20 cells) -0.00215910967243115"
[1] "2022-05-26 00:17:11: Centering and scaling data."
[1] "2022-05-26 00:17:12: Removing genes with no variation."
[1] "2022-05-26 00:17:12: Calculating PCA."
[1] "2022-05-26 00:17:33: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 7.74591785415431"
[1] "5 PCs have eigenvalues larger than 1.5 times null upper bound."
[1] "Storing 10 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

[1] "Mean pseudotime back (~20 cells) 0.00225483200166118"
[1] "Chance of accepted move to equal pseudotime is 0.5"
[1] "Mean pseudotime forward (~20 cells) -0.00225483200166118"
[1] "2022-05-26 00:17:43: Centering and scaling data."
[1] "2022-05-26 00:17:43: Removing genes with no variation."
[1] "2022-05-26 00:17:43: Calculating PCA."
[1] "2022-05-26 00:17:45: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 16.5531277009925"
[1] "4 PCs have eigenvalues larger than 1.5 times null upper bound."
[1] "Storing 8 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

[1] "Mean pseudotime back (~20 cells) 0.00239801845890802"
[1] "Chance of accepted move to equal pseudotime is 0.5"
[1] "Mean pseudotime forward (~20 cells) -0.00239801845890802"
[1] "2022-05-26 00:17:49: Centering and scaling data."
[1] "2022-05-26 00:17:50: Removing genes with no variation."
[1] "2022-05-26 00:17:50: Calculating PCA."
[1] "2022-05-26 00:17:53: Estimating significant PCs."
[1] "Marchenko-Pastur eigenvalue null upper bound: 13.3477783147608"
[1] "5 PCs have eigenvalues larger than 1.5 times null upper bound."
[1] "Storing 10 PCs."

Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

[1] "Mean pseudotime back (~20 cells) 0.00306078733162799"
[1] "Chance of accepted move to equal pseudotime is 0.5"
[1] "Mean pseudotime forward (~20 cells) -0.00306078733162799"

Save

subfix_4_l = list()
for(n in names(sub_urd)){
  subsetfix = urdSubset(sub_urd[[n]], 
                        cells.keep=rownames(sub_urd[[n]]@meta)[rownames(sub_urd[[n]]@meta)%in%rownames(axials_4_l[[n]][["tm"]])])
  root.cells = cellsInCluster(sub_urd[[n]], clustering = "cellclusters", "epen_clus_4")
  
  # Simulate the biased random walks from each tip
  axial.walks = simulateRandomWalksFromTips(subsetfix, tip.group.id="tip.clusters",
                                            root.cells=root.cells, n.per.tip = 25000,
                                            transition.matrix = axials_4_l[[n]][["tm"]], 
                                            root.visits = 1, max.steps = 10000, verbose = F)
  
  # Process the biased random walks into visitation frequencies
  subsetfix = processRandomWalksFromTips(subsetfix, axial.walks, verbose = F)
  
  axials_4_l[[n]][["walks"]] = axial.walks
  subfix_4_l[[n]] = subsetfix
}

Check tip clusters

save(subset_dat_l, floods_4_l, sub_urd, axials_4_l, subfix_4_l,
     file = "data/processed/URD/URD_Div_lists.RData")

Plot trees

max_vals = list()
for(n in names(subfix_4_l)){
  plt1 = plotDim(subfix_4_l[[n]], "cellclusters", plot.title=n)
  plt2 = plotDim(subfix_4_l[[n]], "tip.clusters", plot.title=n)
  print(plt1+plt2)
  
  max_vals[[n]] = sort(tapply(subfix_4_l[[n]]@pseudotime$pseudotime_4,
                    subfix_4_l[[n]]@group.ids$tip.clusters, max))
}
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.
Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> = "none")` instead.

Save trees

cl_tips = list(c("4", "5"), c("4", "2"), c("4", "1"), c("2", "3"), c("1", "2"))
names(cl_tips) = names(subfix_4_l)
tree_plt_l = list()
for(n in names(subfix_4_l)){
  # Load the cells used for each tip into the URD object
  axial.tree = loadTipCells(subfix_4_l[[n]], "tip.clusters")

  # Build the tree
  test = buildTree(axial.tree, pseudotime = "pseudotime_4", tips.use=cl_tips[[n]], 
                  divergence.method = "preference", cells.per.pseudotime.bin = 20,
                  bins.per.pseudotime.window = 8, save.all.breakpoint.info = T,
                  p.thresh=0.001,minimum.visits = 2)
  
  plt1 = plotTree(test, "tip.clusters")
  plt2 = plotTree(test, "cellclusters")
  
  # add labels
  clayout = test@tree$cell.layout
  clayout$ct = test@meta[rownames(clayout),"cellclusters"]
  labdf = data.frame(x = tapply(clayout$x, clayout$ct, mean),
                     y = tapply(clayout$y, clayout$ct, mean))
  labdf$lab = rownames(labdf)
  condx = labdf$y<min(test@tree$segment.pseudotime.limits$end)
  labdf[!condx,"x"] = round(labdf[!condx,"x"], 0)
  labdf[condx,"x"] = 0.5
  
  plt3 = plt2+scale_colour_manual(values = cols_cc, limits = force)+
    ggrepel::geom_label_repel(data = labdf, mapping = aes(x = x, y = y, label = lab), size = 2.2, 
               label.padding = unit(0.1, "lines"), box.padding = unit(0.1, "lines"))+
    labs(title = n)+
    theme(legend.text = element_text(size = 6),
          legend.key.size = unit(0.2, "lines"),
          plot.title = element_text(size = 8))
  print(plt3)
  tree_plt_l[[n]] = plt3
}
[1] "Calculating divergence between 4 and 5 (Pseudotime 0 to 0.348)"
[1] "Joining segments 4 and 5 at pseudotime 0.213 to create segment 6"
[1] "Assigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  4279 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Collapsing short segments."
[1] "Removing singleton segments."
[1] "Reassigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  4279 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Assigning cells to nodes."
[1] "Laying out tree."
[1] "Adding cells to tree."
[1] "Calculating divergence between 4 and 2 (Pseudotime 0 to 0.473)"
Warning in pseudotimeBreakpointByStretch(div.pseudotime, segment.1, segment.2,  :
  No obvious breakpoint between 4 and 2 . Longest stretch of difference is upstream of longest stretch of non-different.
[1] "Joining segments 4 and 2 at pseudotime 0 to create segment 5"
[1] "Assigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  2532 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Collapsing short segments."
[1] "Removing singleton segments."
[1] "Reassigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  2532 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Assigning cells to nodes."
[1] "Laying out tree."
[1] "Adding cells to tree."
[1] "Calculating divergence between 4 and 1 (Pseudotime 0 to 0.522)"
[1] "Joining segments 4 and 1 at pseudotime 0.382 to create segment 5"
[1] "Assigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  2042 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Collapsing short segments."
[1] "Removing singleton segments."
[1] "Reassigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  2042 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Assigning cells to nodes."
[1] "Laying out tree."
[1] "Adding cells to tree."
[1] "Calculating divergence between 2 and 3 (Pseudotime 0 to 0.525)"
Warning in pseudotimeBreakpointByStretch(div.pseudotime, segment.1, segment.2,  :
  No obvious breakpoint between 2 and 3 . Longest stretch of difference is upstream of longest stretch of non-different.
[1] "Joining segments 2 and 3 at pseudotime 0 to create segment 4"
[1] "Assigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  1673 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Collapsing short segments."
[1] "Removing singleton segments."
[1] "Reassigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  1673 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Assigning cells to nodes."
[1] "Laying out tree."
[1] "Adding cells to tree."
[1] "Calculating divergence between 1 and 2 (Pseudotime 0 to 0.613)"
Warning in pseudotimeBreakpointByStretch(div.pseudotime, segment.1, segment.2,  :
  No obvious breakpoint between 1 and 2 . Longest stretch of difference is upstream of longest stretch of non-different.
[1] "Joining segments 1 and 2 at pseudotime 0 to create segment 3"
[1] "Assigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  1310 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Collapsing short segments."
[1] "Removing singleton segments."
[1] "Reassigning cells to segments."
Warning in assignCellsToSegments(object, pseudotime, verbose) :
  1310 cells were not visited by a branch that exists at their pseudotime and were not assigned.
[1] "Assigning cells to nodes."
[1] "Laying out tree."
[1] "Adding cells to tree."

Load RNA velocity pseudotime

for(n in names(tree_plt_l)){
  pdf(paste0("results/pseudotime/Div_", n, "_tree.pdf"), height = 3.5, width = 4)
  print(tree_plt_l[[n]])
  dev.off()
}

Plot pseudotime comparisons

glut_dat_df = read.csv("results/RNAvelocity_Div/heatmap_gen/div_dat_df.csv", 
                       header = T, row.names = 1)
newcellnames = rownames(glut_dat_df)
newcellnames = gsub("_2_wpi", "-1_1", newcellnames)
newcellnames = gsub("_4_wpi", "-1_2", newcellnames)
newcellnames = gsub("_6_wpi", "-1_3", newcellnames)
newcellnames = gsub("_8_wpi", "-1_4", newcellnames)
newcellnames = gsub("_12_wpi", "-1_5", newcellnames)
newcellnames = gsub("_1_wpi_pos", "-1", newcellnames)
rownames(glut_dat_df) = newcellnames

Save comparisons

plt_comp_l = list()
for(n in names(subfix_l)){
  plot_df = merge(glut_dat_df, subfix_4_l[[n]]@pseudotime, by = 0)
  cpe = round(cor(plot_df$newpt, plot_df$pseudotime_4, method = "pe"), 2)
  
  plt_comp_l[[n]] = ggplot(plot_df, aes(x = newpt, y = pseudotime_4, 
                                        colour = pred_ctall))+
    geom_point()+
    scale_colour_manual(values = cols_cc, limits = force)+
    #scale_size_manual(values = c(0.5, 0.95, 1.35))+
    labs(title = n, subtitle = paste0("PCC = ", cpe), 
         x = "Pseudotime (scVelo)", y = "Pseudotime (URD)")+
    theme_classic()+
    theme(aspect.ratio = 1,
          axis.text = element_text(colour = "black", size = 6),
          axis.title = element_text(size = 7),
          plot.title = element_text(size = 8),
          plot.subtitle = element_text(size = 7),
          legend.text = element_text(size = 6),
          legend.title = element_text(size = 7),
          legend.key.size = unit(0.35, "cm"))
}

Make proportion plots

for(n in names(plt_comp_l)){
  pdf(paste0("results/pseudotime/Div_", n, "_comp.pdf"), height = 4, width = 4.5)
  print(plt_comp_l[[n]])
  dev.off()
}

Save proportions

smo_prop_list = list()
for(n in names(subfix_4_l)){
  # subset data
  submeta = data.frame("pseudotime" = subfix_4_l[[n]]@pseudotime$pseudotime_4,
                       "cellclusters" = subfix_4_l[[n]]@meta$cellclusters)
  
  med_pt_cc = sort(tapply(submeta$pseudotime, submeta$cellclusters, median))
  
  lt_bins = cut(submeta$pseudotime, 100) # 100 equally-sized bins
  plot_df = data.frame(bins = lt_bins, 
                       cst = as.character(submeta$cellclusters))
  tab_df = table(plot_df$bins, plot_df$cst)
  
  # remove cell types that are too rare (<5%)
  tab_df = reshape2::melt(tab_df/rowSums(tab_df))
  usecl = tapply(tab_df$value, tab_df$Var2, function(x) any(x>0.05))
  plot_df = plot_df[plot_df$cst %in% names(usecl)[usecl],]
  tab_df = table(plot_df$bins,plot_df$cst)
  
  # normalise by cell type abundance
  prop_w = prop.table(table(plot_df$cst))
  tab_df = t(apply(tab_df, 1, function(x) x/prop_w[colnames(tab_df)]))
  
  # reshape
  tab_df = reshape2::melt(tab_df/rowSums(tab_df))
  tab_df$Var2 = as.character(tab_df$Var2)
  
  # prevent discontinuity by copying the previous column (likely not happening)
  tab_df = tab_df[order(tab_df$Var1, decreasing = F),]
  for(i in unique(tab_df$Var1)){
    if(any(is.nan(tab_df$value[tab_df$Var1==i]))){
      tab_df$value[tab_df$Var1==i] = prev
    }
    prev = tab_df$value[tab_df$Var1==i]
  }
  
  # smoothen the proportions (and force constrain to 0-1)
  tab_df2 = tab_df
  tab_df2$value2 = tab_df2$value
  for(i in unique(tab_df2$Var2)){
    fff = loess(value~as.numeric(Var1), data = tab_df2[tab_df2$Var2==i,], 
                span = 0.5)
    pred = predict(fff)
    pred[pred>1] = 1
    pred[pred<0] = 0
    tab_df2$value2[tab_df2$Var2==i] = pred
  }
  
  # force constrain each interval to 0-1 by doing proportion
  for(i in unique(tab_df2$Var1)){
    tab_df2$value2[tab_df2$Var1==i] = tab_df2$value2[tab_df2$Var1==i]/sum(tab_df2$value2[tab_df2$Var1==i])
  }
  
  tab_df2$major = unlist(lapply(strsplit(tab_df2$Var2, "_"), function(x) x[1]))
  res = list()
  for(nnn in unique(tab_df2$Var1)){
    ss=tapply(tab_df2[tab_df2$Var1==nnn,"value2"], tab_df2[tab_df2$Var1==nnn,"major"], sum)
    res[[nnn]] = which.max(ss)
  }

  smo_prop_list[[n]] = ggplot(tab_df2, aes(x = Var1, y = value2, group = Var2, fill = Var2))+
    geom_area()+
    scale_y_continuous(expand = c(0,0))+
    scale_fill_manual(values = cols_cc[names(cols_cc) %in% tab_df2$Var2])+
    labs(x = "Bins", y = "Proportion", fill = "Cell type")+
    theme_classic()+
    theme(axis.text.x = element_blank(),
          axis.text.y = element_text(size = 6.5, colour = "black"),
          axis.ticks.x = element_blank(),
          axis.line = element_blank(),
          axis.title = element_text(size = 7),
          legend.text = element_text(size = 6),
          legend.title = element_text(size = 7),
          legend.key.size = unit(0.4, "cm"))
}
LS0tCnRpdGxlOiAiUHNldWRvdGltZSBmb3IgRGl2LXNlcSIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKCgojIEdlbmVyYWwgc2V0dXAKU2V0dXAgY2h1bmsKCmBgYHtyLCBzZXR1cCwgaW5jbHVkZT1GQUxTRX0Ka25pdHI6Om9wdHNfY2h1bmskc2V0KGZpZy53aWR0aCA9IDgpCmtuaXRyOjpvcHRzX2tuaXQkc2V0KHJvb3QuZGlyID0gbm9ybWFsaXplUGF0aCgiLi4iKSkKa25pdHI6Om9wdHNfa25pdCRnZXQoInJvb3QuZGlyIikKYGBgCgpMb2FkIGxpYnJhcmllcwoKYGBge3J9CmxpYnJhcnkoU2V1cmF0KQpsaWJyYXJ5KGhhcm1vbnkpCmxpYnJhcnkoZ2dwbG90MikKbGlicmFyeShNYXRyaXgpCmxpYnJhcnkobWdjdikKbGlicmFyeShmb3JlYWNoKQpsaWJyYXJ5KGRvUGFyYWxsZWwpCmxpYnJhcnkocGFyYWxsZWwpCmxpYnJhcnkoVVJEKQpsaWJyYXJ5KHNsaW5nc2hvdCkKbGlicmFyeShSQ29sb3JCcmV3ZXIpCmBgYAoKU2V0IGNvbG91cnMgZm9yIGNlbGwgdHlwZXMgYW5kIHJlZ2lvbnMKCmBgYHtyfQptZXRhID0gcmVhZC5jc3YoImRhdGEvYW5ub3RhdGlvbnMvYXhvbG90bF9hbGxfdW1ldGEuY3N2IiwgCiAgICAgICAgICAgICAgICBoZWFkZXIgPSBULCByb3cubmFtZXMgPSAxKQpjb2xzX2NjID0gYygKI2VwZW4KIiMxMjQwMGMiLCAiIzJkNjYyNCIsIiMxZDRmMTUiLCAiIzE3NDcxMSIsICIjMmQ2NjI0IiwgIiMzZDdmMzMiLCAiIzNiN2IzMCIsICIjNDY4YjNiIiwgIiM0Zjk4NDMiLCIjNWRhZTUwIiwgIiM2NmJiNTgiLCAiIzcyY2Q2NCIsICIjMzA2YTI2IiwgIiM3OGQ2NjkiLCAiIzgxZTQ3MiIsCiNnYWJhCiIjNzAwMjA5IiwgIiM3NTA5MGUiLCIjN2EwZjEzIiwgIiM4MDE1MTciLCAiIzg1MWExYiIsICIjOGExZjFmIiwgIiM5MDI0MjMiLCAiIzk1MjkyNyIsICIjOWEyZDJjIiwiI2EwMzIzMCIsICIjYTUzNjM0IiwgIiNhYTNhMzkiLCAiI2IwM2YzZCIsIiNiNTQzNDIiLCAiI2JhNDg0NiIsICIjYzA0YzRiIiwgIiNjNTUwNGYiLCAiI2NhNTU1NCIsICIjZDA1OTU5IiwgIiNkNTVlNWUiLCIjNzMwNTBjIiwgIiM3ODBjMTEiLCIjOGQyMjIxIiwgIiM5ODJiMmEiLCIjYTIzNDMyIiwgIiNhODM4MzciLCAiI2IyNDEzZiIsICIjYjg0NTQ0IiwgIiNiZDRhNDkiLCAiI2M4NTM1MiIsICMiI2NkNTc1NiIsCiNnbHV0CiIjMDU0Njc0IiwgIiMxMzRkN2IiLCIjMWQ1NDgxIiwgIiMyNjVhODgiLCAiIzJlNjE4ZSIsICIjNzNhNGNiIiwgIiMzNjY5OTUiLCAiIzNlNzA5YyIsICIjNDY3N2EyIiwiIzRkN2VhOSIsICIjNTU4NmIwIiwgIiM1YzhkYjciLCAiIzY0OTViZCIsIiM2YjljYzQiLCAiIzdiYWNkMiIsICIjOGViZmU0IiwgIiM5NmM3ZWIiLCAiIzllY2ZmMiIsICIjMTg1MDdlIiwgIiMxODUwN2UiLCIjMmE1ZThiIiwgIiM0OTdiYTYiLCIjNTg4OWIzIiwgIiM2ZmEwYzgiLCIjN2ZhZmQ2IiwgIiM2MDkxYmEiLCAiIzUxODJhYyIsICIjM2E2Yzk4IiwgIiNhNmQ3ZjkiLAojbnBjCiIjZmZiMTIwIiwgIiNmZWI3MmEiLCIjZmRiYzM0IiwgIiNmY2MxM2QiLCAiI2ZiYzc0NSIsICIjZmFjYzRlIiwgIiNmOWQxNTYiLCAiI2Y4ZDY1ZiIsICIjZjhkYTY4IiwiI2Y3ZGY3MCIsICIjZjdlNDc5IiwgIiNmN2U4ODIiLCAiI2Y3ZWQ4YSIsICIjZjdmMTkzIiwgIiNlY2EzMTkiCikKY2NuYW1lcyA9IHVuaXF1ZShzb3J0KG1ldGEkY2VsbGNsdXN0ZXJzKSkKbmFtZXMoY29sc19jYykgPSBjKGNjbmFtZXNbZ3JlcGwoImVwZW4iLCBjY25hbWVzKV0sIGNjbmFtZXNbZ3JlcGwoIkdBQkEiLCBjY25hbWVzKV0sY2NuYW1lc1tncmVwbCgiZ2x1dCIsIGNjbmFtZXMpXSxjY25hbWVzW2dyZXBsKCJucGMiLCBjY25hbWVzKV0pCgpyZWdfY29scyA9IGMoIm90aGVyL3Vua25vd25fcHJlZCIgPSAiI0M3Q0NDNyIsIAogICAgICAgICAgICAgIm1lZGlhbCIgPSAiIzUyMTY4RCIsICJtZWRpYWxfcHJlZCIgPSAiIzY2MUNCMCIsIAogICAgICAgICAgICAgImRvcnNhbCIgPSAiI0M1NjAwNyIsICJkb3JzYWxfcHJlZCIgPSAiI0VENzMwNyIsIAogICAgICAgICAgICAgImxhdGVyYWwiID0gIiMxMTgzOTIiLCAibGF0ZXJhbF9wcmVkIiA9ICIjMTZBM0I2IikKcmVnX2NvbHNfc2ltcCA9IGMoIm1lZGlhbCIgPSAiIzUyMTY4RCIsICJkb3JzYWwiID0gIiNDNTYwMDciLCAibGF0ZXJhbCIgPSAiIzExODM5MiIpCgpjb2xvcnMgPC0gY29sb3JSYW1wUGFsZXR0ZShicmV3ZXIucGFsKDExLCdTcGVjdHJhbCcpWy02XSkoMTAwKQpgYGAKCgoKTG9hZCBkYXRhCgpgYGB7cn0KZGl2X3NyYXQgPSByZWFkUkRTKCJkYXRhL2V4cHJlc3Npb24vYXhvbG90bF9yZWNsdXN0L0VkdV8xXzJfNF82XzhfMTJfZmlsX2hpZ2h2YXJmZWF0LlJEUyIpCnByZWRfbWV0YSA9IHJlYWQuY3N2KCJyZXN1bHRzL0Rpdi1zZXEvZGl2c2VxX3ByZWRpY3RlZF9tZXRhZGF0YS5jc3YiLCAKICAgICAgICAgICAgICAgICAgICAgaGVhZGVyID0gVCwgcm93Lm5hbWVzID0gMSlbLDQ6NV0KY29sbmFtZXMocHJlZF9tZXRhKSA9IGMoInJlZ2lvbiIsICJjZWxsY2x1c3RlcnMiKQpkaXZfc3JhdCA9IEFkZE1ldGFEYXRhKGRpdl9zcmF0LCBtZXRhZGF0YSA9IHByZWRfbWV0YSkKYGBgCgoKCiMgQ2FsY3VsYXRlIHBzZXVkb3RpbWUKIyMgUHJvY2VzcyBkYXRhCkRlZmluZSBncm91cHMKCmBgYHtyfQpjb21tb25fY2VsbHMgPSBjKCJlcGVuX2NsdXNfMyIsICJlcGVuX2NsdXNfNCIpCnJfY2VsbHMgPSBsaXN0KCJoaXBwb2NhbXB1cyIgPSBjKCJnbHV0X1NVQlNFVF8wIiwgImdsdXRfU1VCU0VUXzQiLCAiZ2x1dF9TVUJTRVRfNSIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJnbHV0X1NVQlNFVF83IiwgIm5wY19TVUJTRVRfNyIsIm5wY19TVUJTRVRfMTEiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAiZ2x1dF9TVUJTRVRfMTQiLCAiZ2x1dF9TVUJTRVRfMTUiLCAiZ2x1dF9TVUJTRVRfMTIiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAiZ2x1dF9TVUJTRVRfMTYiLCAiZ2x1dF9TVUJTRVRfMTciLCAiZ2x1dF9TVUJTRVRfMyIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJucGNfU1VCU0VUXzEiLCAibnBjX1NVQlNFVF8zIiwgIm5wY19TVUJTRVRfOSIpLAogICAgImxjIiA9IGMoImdsdXRfU1VCU0VUXzkiLCJucGNfU1VCU0VUXzciLCJucGNfU1VCU0VUXzQiLCAiZ2x1dF9TVUJTRVRfMjEiLCAiZ2x1dF9TVUJTRVRfMiIsCiAgICAgICAgICAgICAiZ2x1dF9TVUJTRVRfMjUiLCJucGNfU1VCU0VUXzE0IiwibnBjX1NVQlNFVF8wIiksCiAgICAiY2wxMTEiID0gYygiZ2x1dF9TVUJTRVRfMSIsICJnbHV0X1NVQlNFVF8xMSIsICJucGNfU1VCU0VUXzIiLCJucGNfU1VCU0VUXzQiLAogICAgICAgICAgICAgICAgIm5wY19TVUJTRVRfNyIsICJucGNfU1VCU0VUXzEzIiksCiAgICAiZW9tZXMiID0gYygibnBjX1NVQlNFVF83IiwibnBjX1NVQlNFVF80IiwgImdsdXRfU1VCU0VUXzEwIiwiZ2x1dF9TVUJTRVRfMjIiKSwKICAgICJjbDg2MjBfZXAiID0gYygiZ2x1dF9TVUJTRVRfOCIsImdsdXRfU1VCU0VUXzYiLCJnbHV0X1NVQlNFVF8yMCIsICJlcGVuX2NsdXNfMSIsCiAgICAgICAgICAgICAgICAgICAgImVwZW5fY2x1c183IiwgImVwZW5fY2x1c18xNCIpKQpgYGAKClN1YnNldCBkYXRhCgpgYGB7cn0Kc3ViX3NyYXQgPSBsaXN0KCkKZm9yKG4gaW4gbmFtZXMocl9jZWxscykpewogIHN1Yl9zcmF0W1tuXV0gPSBkaXZfc3JhdFssZGl2X3NyYXQkY2VsbGNsdXN0ZXJzICVpbiUgYyhjb21tb25fY2VsbHMsIHJfY2VsbHNbW25dXSldCn0KYGBgCgpRdWljayBwcm9jZXNzaW5nCgpgYGB7cn0KZm9yKG4gaW4gbmFtZXMoc3ViX3NyYXQpKXsKICBwcmludChuKQogIHN1Yl9zcmF0W1tuXV0gPSBOb3JtYWxpemVEYXRhKHN1Yl9zcmF0W1tuXV0pCiAgc3ViX3NyYXRbW25dXSA9IEZpbmRWYXJpYWJsZUZlYXR1cmVzKHN1Yl9zcmF0W1tuXV0pCiAgc3ViX3NyYXRbW25dXSA9IFNjYWxlRGF0YShzdWJfc3JhdFtbbl1dLCB2YXJzLnRvLnJlZ3Jlc3MgPSBjKCJuQ291bnRfUk5BIiksIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgZmVhdHVyZXMgPSBWYXJpYWJsZUZlYXR1cmVzKG9iamVjdCA9IHN1Yl9zcmF0W1tuXV0pKQp9CmBgYAoKCgojIFVSRApQcmVwYXJlIFVSRCBvYmplY3RzCgpgYGB7cn0Kc3ViX3VyZCA9IGxpc3QoKQpmb3IobiBpbiBuYW1lcyhzdWJfc3JhdCkpewogIHN1Yl91cmRbW25dXSA9IHN1Yl9zcmF0W1tuXV0KICBzdWJfdXJkW1tuXV0gPSBjcmVhdGVVUkQoY291bnQuZGF0YSA9IHN1Yl91cmRbW25dXUBhc3NheXMkUk5BQGNvdW50cywgCiAgICAgICAgICAgICAgICAgICAgICAgICAgIG1ldGEgPSBzdWJfdXJkW1tuXV1AbWV0YS5kYXRhLCBtaW4uY2VsbHM9MywgbWluLmNvdW50cz0zKQogIHN1Yl91cmRbW25dXUBncm91cC5pZHMkY2VsbGNsdXN0ZXJzID0gc3ViX3VyZFtbbl1dQG1ldGEkY2VsbGNsdXN0ZXJzCiAgc3ViX3VyZFtbbl1dQGdyb3VwLmlkcyRzdWJjbGFzc2VzID0gc3ViX3VyZFtbbl1dQG1ldGEkc3ViY2xhc3Nlcwp9CmBgYAoKQ2FsY3VsYXRlIERScwoKYGBge3J9CmZvcihuIGluIG5hbWVzKHN1Yl91cmQpKXsKICBzdWJfdXJkW1tuXV1AdmFyLmdlbmVzID0gZmluZFZhcmlhYmxlR2VuZXMoc3ViX3VyZFtbbl1dLCBzZXQub2JqZWN0LnZhci5nZW5lcz1GLGRvLnBsb3Q9RiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGRpZmZDVi5jdXRvZmY9MC4zLCBtZWFuLm1pbj0uMDA1LCBtZWFuLm1heD0xMDApCgogIHN1Yl91cmRbW25dXSA9IGNhbGNQQ0Eoc3ViX3VyZFtbbl1dLCBtcC5mYWN0b3IgPSAyKQogIHBjU0RQbG90KHN1Yl91cmRbW25dXSkKICBzdWJfdXJkW1tuXV0gPSBjYWxjVHNuZShvYmplY3QgPSBzdWJfdXJkW1tuXV0pCiAgcGxvdERpbShzdWJfdXJkW1tuXV0sICJjZWxsY2x1c3RlcnMiKQogIAogIHN1Yl91cmRbW25dXSA9IGNhbGNETShzdWJfdXJkW1tuXV0sIGtubiA9IDE1MCwgc2lnbWE9MTYpCiAgCiAgcGxvdERpbUFycmF5KHN1Yl91cmRbW25dXSwgcmVkdWN0aW9uLnVzZSA9ICJkbSIsIGRpbXMudG8ucGxvdCA9IDE6OCwgCiAgICAgICAgICAgICBvdXRlci50aXRsZSA9ICJEaWZmdXNpb24gTWFwIChTaWdtYSAxNiwgMTUwIE5Ocyk6IFN0YWdlIiwgbGFiZWw9ImNlbGxjbHVzdGVycyIsCiAgICAgICAgICAgICBwbG90LnRpdGxlPSIiLCBsZWdlbmQ9RikKCiAgcGxvdERpbShzdWJfdXJkW1tuXV0sICJjZWxsY2x1c3RlcnMiLCB0cmFuc2l0aW9ucy5wbG90ID0gMTAwMDApCn0KYGBgCgpQc2V1ZG90aW1lIGNhbGN1bGF0aW9uCgpgYGB7cn0KZW5kX2NjID0gbGlzdCgiaGlwcG9jYW1wdXMiID0gYygiZ2x1dF9TVUJTRVRfMCIsICJnbHV0X1NVQlNFVF83IiksCiAgICAgICAgICAgICAgImxjIiA9IGMoImdsdXRfU1VCU0VUXzkiLCJnbHV0X1NVQlNFVF8yIiksCiAgICAgICAgICAgICAgImNsMTExIiA9IGMoImdsdXRfU1VCU0VUXzEiLCAiZ2x1dF9TVUJTRVRfMTEiKSwKICAgICAgICAgICAgICAiZW9tZXMiID0gYygiZ2x1dF9TVUJTRVRfMTAiLCJnbHV0X1NVQlNFVF8yMiIpLAogICAgICAgICAgICAgICJjbDg2MjBfZXAiID0gYygiZ2x1dF9TVUJTRVRfOCIsImdsdXRfU1VCU0VUXzYiKSkKc3Vic2V0X2RhdF9sID0gbGlzdCgpCmZsb29kc180X2wgPSBsaXN0KCkKZm9yKG4gaW4gbmFtZXMoc3ViX3VyZCkpewogIHByaW50KG4pCiAgIyBIZXJlIHdlIHVzZSBhbGwgY2VsbHMgZnJvbSB0aGUgZmlyc3Qgc3RhZ2UgYXMgdGhlIHJvb3QKICByb290LmNlbGxzID0gY2VsbHNJbkNsdXN0ZXIoc3ViX3VyZFtbbl1dLCBjbHVzdGVyaW5nID0gImNlbGxjbHVzdGVycyIsICJlcGVuX2NsdXNfNCIpCiAgCiAgIyBUaGVuIHdlIHJ1biAnZmxvb2QnIHNpbXVsYXRpb25zCiAgZmxvb2RzID0gZmxvb2RQc2V1ZG90aW1lKHN1Yl91cmRbW25dXSwgcm9vdC5jZWxscyA9IHJvb3QuY2VsbHMsIG49MjUwLAogICAgICAgICAgICAgICAgICAgICAgICAgICBtaW5pbXVtLmNlbGxzLmZsb29kZWQgPSAyLCB2ZXJib3NlPUYpCiAgIyBmaXggZXhjZXNzaXZlIE5BCiAgZmxvb2RzID0gZmxvb2RzWyxjb2xTdW1zKCFpcy5uYShmbG9vZHMpKT4xMDAwXQogIAogICMgVGhlIHdlIHByb2Nlc3MgdGhlIHNpbXVsYXRpb25zIGludG8gYSBwc2V1ZG90aW1lCiAgc3ViX3VyZFtbbl1dID0gZmxvb2RQc2V1ZG90aW1lUHJvY2VzcyhzdWJfdXJkW1tuXV0sIGZsb29kcywgZmxvb2RzLm5hbWU9InBzZXVkb3RpbWVfNCIpCiAgCiAgcHNldWRvdGltZVBsb3RTdGFiaWxpdHlPdmVyYWxsKHN1Yl91cmRbW25dXSkKICAKICBwbG90RGlzdHMoc3ViX3VyZFtbbl1dLCAicHNldWRvdGltZV80IiwgImNlbGxjbHVzdGVycyIsIHBsb3QudGl0bGU9IlBzZXVkb3RpbWUgYnkgc3RhZ2UiKQogIAogICMgQ3JlYXRlIGEgc3Vic2V0dGVkIG9iamVjdCBvZiBqdXN0IHRob3NlIGNlbGxzIGZyb20gdGhlIGZpbmFsIHN0YWdlCiAgc3Vic2V0ZGF0ID0gdXJkU3Vic2V0KHN1Yl91cmRbW25dXSwgCiAgICAgICAgICAgICAgICAgICAgICAgIGNlbGxzLmtlZXA9Y2VsbHNJbkNsdXN0ZXIoc3ViX3VyZFtbbl1dLCAiY2VsbGNsdXN0ZXJzIiwgZW5kX2NjW1tuXV0pKQogIHN1YnNldF9kYXRfbFtbbl1dID0gc3Vic2V0ZGF0CiAgZmxvb2RzXzRfbFtbbl1dID0gZmxvb2RzCn0KYGBgCgpTYXZlCgpgYGB7cn0Kc2F2ZShzdWJzZXRfZGF0X2wsIHN1Yl91cmQsIGZsb29kc180X2wsIGZpbGUgPSAiZGF0YS9wcm9jZXNzZWQvVVJEL1VSRF9EaXZfbGlzdHMuUkRhdGEiKQpgYGAKCkRldGVybWluZSBuZXcgaWRlbnRpdGllcyBmb3IgZW5kIHRpcHMKCmBgYHtyfQpheGlhbHNfNF9sID0gbGlzdCgpCmZvcihuIGluIG5hbWVzKHN1Yl91cmQpKXsKICAjIFVzZSB0aGUgdmFyaWFibGUgZ2VuZXMgdGhhdCB3ZXJlIGNhbGN1bGF0ZWQgb25seSBvbiB0aGUgZmluYWwgZ3JvdXAgb2Ygc3RhZ2VzICh3aGljaAogICMgY29udGFpbiB0aGUgbGFzdCBzdGFnZSkuCiAgc3Vic2V0X2RhdF9sW1tuXV1AdmFyLmdlbmVzID0gc3ViX3VyZFtbbl1dQHZhci5nZW5lc1tzdWJfdXJkW1tuXV1AdmFyLmdlbmVzICVpbiUgcm93bmFtZXMoc3Vic2V0X2RhdF9sW1tuXV1AbG9ndXB4LmRhdGEpXQogIAogICMgQ2FsY3VsYXRlIFBDQSBhbmQgdFNORQogIHN1YnNldF9kYXRfbFtbbl1dID0gY2FsY1BDQShzdWJzZXRfZGF0X2xbW25dXSwgbXAuZmFjdG9yID0gMS41KQogIHBjU0RQbG90KHN1YnNldF9kYXRfbFtbbl1dKQogIAogIHNldC5zZWVkKDIwKQogIHN1YnNldGRhdCA9IGNhbGNUc25lKHN1YnNldF9kYXRfbFtbbl1dKQogIAogICMgQ2FsY3VsYXRlIGdyYXBoIGNsdXN0ZXJpbmcgb2YgdGhlc2UgY2VsbHMKICBzdWJzZXRfZGF0X2xbW25dXSA9IGdyYXBoQ2x1c3RlcmluZyhzdWJzZXRfZGF0X2xbW25dXSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbnVtLm5uID0gNTAsIGRvLmphY2NhcmQ9VCwgbWV0aG9kPSJMb3V2YWluIikKICBwbG90RGltKHN1YnNldF9kYXRfbFtbbl1dLCAiTG91dmFpbi01MCIsIAogICAgICAgICAgcGxvdC50aXRsZSA9ICJMb3V2YWluICg1MCBOTikgZ3JhcGggY2x1c3RlcmluZyIsIHBvaW50LnNpemU9MykKICBwbG90RGltKHN1YnNldF9kYXRfbFtbbl1dLCAiY2VsbGNsdXN0ZXJzIiwgcGxvdC50aXRsZSA9ICJjZWxsY2x1c3RlcnMiLCBwb2ludC5zaXplPTMpCiAgCiAgCiAgIyBDb3B5IGNsdXN0ZXIgaWRlbnRpdGllcyBmcm9tIGF4aWFsLjZzb21pdGUgb2JqZWN0IHRvIGEgbmV3IGNsdXN0ZXJpbmcgKCJ0aXAuY2x1c3RlcnMiKSBpbiB0aGUgZnVsbCBheGlhbCBvYmplY3QuCiAgc3ViX3VyZFtbbl1dQGdyb3VwLmlkc1tyb3duYW1lcyhzdWJzZXRfZGF0X2xbW25dXUBncm91cC5pZHMpLCAidGlwLmNsdXN0ZXJzIl0gPSBzdWJzZXRfZGF0X2xbW25dXUBncm91cC5pZHMkYExvdXZhaW4tNTBgCiAgCiAgIyBEZXRlcm1pbmUgdGhlIHBhcmFtZXRlcnMgb2YgdGhlIGxvZ2lzdGljIHVzZWQgdG8gYmlhcyB0aGUgdHJhbnNpdGlvbiBwcm9iYWJpbGl0aWVzLiAKICAjIFRoZSBwcm9jZWR1cmUgaXMgcmVsYXRpdmVseSByb2J1c3QgdG8gdGhpcyBwYXJhbWV0ZXIsIGJ1dCB0aGUgY2VsbCBudW1iZXJzIG1heSBuZWVkIHRvIGJlCiAgIyBtb2RpZmllZCBmb3IgbGFyZ2VyIG9yIHNtYWxsZXIgZGF0YSBzZXRzLgogIGF4aWFsLnB0bG9naXN0aWMgPSBwc2V1ZG90aW1lRGV0ZXJtaW5lTG9naXN0aWMoc3ViX3VyZFtbbl1dLCAicHNldWRvdGltZV80IiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIG9wdGltYWwuY2VsbHMuZm9yd2FyZD0yMCwgbWF4LmNlbGxzLmJhY2s9MjAsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBkby5wbG90ID0gVCkKICAKICBheGlhbC5iaWFzZWQudG0gPSBhcy5tYXRyaXgocHNldWRvdGltZVdlaWdodFRyYW5zaXRpb25NYXRyaXgoc3ViX3VyZFtbbl1dLCAicHNldWRvdGltZV80IiwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgbG9naXN0aWMucGFyYW1zPWF4aWFsLnB0bG9naXN0aWMpKQogIAogIGF4aWFsc180X2xbW25dXSA9IGxpc3QoImxvZyIgPSBheGlhbC5wdGxvZ2lzdGljLCAidG0iID0gYXhpYWwuYmlhc2VkLnRtKQp9CmBgYAoKUmFuZG9tIHdhbGsgc2ltdWxhdGlvbnMKCmBgYHtyfQpzdWJmaXhfNF9sID0gbGlzdCgpCmZvcihuIGluIG5hbWVzKHN1Yl91cmQpKXsKICBzdWJzZXRmaXggPSB1cmRTdWJzZXQoc3ViX3VyZFtbbl1dLCAKICAgICAgICAgICAgICAgICAgICAgICAgY2VsbHMua2VlcD1yb3duYW1lcyhzdWJfdXJkW1tuXV1AbWV0YSlbcm93bmFtZXMoc3ViX3VyZFtbbl1dQG1ldGEpJWluJXJvd25hbWVzKGF4aWFsc180X2xbW25dXVtbInRtIl1dKV0pCiAgcm9vdC5jZWxscyA9IGNlbGxzSW5DbHVzdGVyKHN1Yl91cmRbW25dXSwgY2x1c3RlcmluZyA9ICJjZWxsY2x1c3RlcnMiLCAiZXBlbl9jbHVzXzQiKQogIAogICMgU2ltdWxhdGUgdGhlIGJpYXNlZCByYW5kb20gd2Fsa3MgZnJvbSBlYWNoIHRpcAogIGF4aWFsLndhbGtzID0gc2ltdWxhdGVSYW5kb21XYWxrc0Zyb21UaXBzKHN1YnNldGZpeCwgdGlwLmdyb3VwLmlkPSJ0aXAuY2x1c3RlcnMiLAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHJvb3QuY2VsbHM9cm9vdC5jZWxscywgbi5wZXIudGlwID0gMjUwMDAsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgdHJhbnNpdGlvbi5tYXRyaXggPSBheGlhbHNfNF9sW1tuXV1bWyJ0bSJdXSwgCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgcm9vdC52aXNpdHMgPSAxLCBtYXguc3RlcHMgPSAxMDAwMCwgdmVyYm9zZSA9IEYpCiAgCiAgIyBQcm9jZXNzIHRoZSBiaWFzZWQgcmFuZG9tIHdhbGtzIGludG8gdmlzaXRhdGlvbiBmcmVxdWVuY2llcwogIHN1YnNldGZpeCA9IHByb2Nlc3NSYW5kb21XYWxrc0Zyb21UaXBzKHN1YnNldGZpeCwgYXhpYWwud2Fsa3MsIHZlcmJvc2UgPSBGKQogIAogIGF4aWFsc180X2xbW25dXVtbIndhbGtzIl1dID0gYXhpYWwud2Fsa3MKICBzdWJmaXhfNF9sW1tuXV0gPSBzdWJzZXRmaXgKfQpgYGAKClNhdmUKCmBgYHtyfQpzYXZlKHN1YnNldF9kYXRfbCwgZmxvb2RzXzRfbCwgc3ViX3VyZCwgYXhpYWxzXzRfbCwgc3ViZml4XzRfbCwKICAgICBmaWxlID0gImRhdGEvcHJvY2Vzc2VkL1VSRC9VUkRfRGl2X2xpc3RzLlJEYXRhIikKYGBgCgpDaGVjayB0aXAgY2x1c3RlcnMKCmBgYHtyLCBmaWcuaGVpZ2h0PTMuNSwgZmlnLndpZHRoPTh9Cm1heF92YWxzID0gbGlzdCgpCmZvcihuIGluIG5hbWVzKHN1YmZpeF80X2wpKXsKICBwbHQxID0gcGxvdERpbShzdWJmaXhfNF9sW1tuXV0sICJjZWxsY2x1c3RlcnMiLCBwbG90LnRpdGxlPW4pCiAgcGx0MiA9IHBsb3REaW0oc3ViZml4XzRfbFtbbl1dLCAidGlwLmNsdXN0ZXJzIiwgcGxvdC50aXRsZT1uKQogIHByaW50KHBsdDErcGx0MikKICAKICBtYXhfdmFsc1tbbl1dID0gc29ydCh0YXBwbHkoc3ViZml4XzRfbFtbbl1dQHBzZXVkb3RpbWUkcHNldWRvdGltZV80LAogICAgICAgICAgICAgICAgICAgIHN1YmZpeF80X2xbW25dXUBncm91cC5pZHMkdGlwLmNsdXN0ZXJzLCBtYXgpKQp9CmBgYAoKUGxvdCB0cmVlcwoKYGBge3IsIGZpZy5oZWlnaHQ9My41LCBmaWcud2lkdGg9NH0KY2xfdGlwcyA9IGxpc3QoYygiNCIsICI1IiksIGMoIjQiLCAiMiIpLCBjKCI0IiwgIjEiKSwgYygiMiIsICIzIiksIGMoIjEiLCAiMiIpKQpuYW1lcyhjbF90aXBzKSA9IG5hbWVzKHN1YmZpeF80X2wpCnRyZWVfcGx0X2wgPSBsaXN0KCkKZm9yKG4gaW4gbmFtZXMoc3ViZml4XzRfbCkpewogICMgTG9hZCB0aGUgY2VsbHMgdXNlZCBmb3IgZWFjaCB0aXAgaW50byB0aGUgVVJEIG9iamVjdAogIGF4aWFsLnRyZWUgPSBsb2FkVGlwQ2VsbHMoc3ViZml4XzRfbFtbbl1dLCAidGlwLmNsdXN0ZXJzIikKCiAgIyBCdWlsZCB0aGUgdHJlZQogIHRlc3QgPSBidWlsZFRyZWUoYXhpYWwudHJlZSwgcHNldWRvdGltZSA9ICJwc2V1ZG90aW1lXzQiLCB0aXBzLnVzZT1jbF90aXBzW1tuXV0sIAogICAgICAgICAgICAgICAgICBkaXZlcmdlbmNlLm1ldGhvZCA9ICJwcmVmZXJlbmNlIiwgY2VsbHMucGVyLnBzZXVkb3RpbWUuYmluID0gMjAsCiAgICAgICAgICAgICAgICAgIGJpbnMucGVyLnBzZXVkb3RpbWUud2luZG93ID0gOCwgc2F2ZS5hbGwuYnJlYWtwb2ludC5pbmZvID0gVCwKICAgICAgICAgICAgICAgICAgcC50aHJlc2g9MC4wMDEsbWluaW11bS52aXNpdHMgPSAyKQogIAogIHBsdDEgPSBwbG90VHJlZSh0ZXN0LCAidGlwLmNsdXN0ZXJzIikKICBwbHQyID0gcGxvdFRyZWUodGVzdCwgImNlbGxjbHVzdGVycyIpCiAgCiAgIyBhZGQgbGFiZWxzCiAgY2xheW91dCA9IHRlc3RAdHJlZSRjZWxsLmxheW91dAogIGNsYXlvdXQkY3QgPSB0ZXN0QG1ldGFbcm93bmFtZXMoY2xheW91dCksImNlbGxjbHVzdGVycyJdCiAgbGFiZGYgPSBkYXRhLmZyYW1lKHggPSB0YXBwbHkoY2xheW91dCR4LCBjbGF5b3V0JGN0LCBtZWFuKSwKICAgICAgICAgICAgICAgICAgICAgeSA9IHRhcHBseShjbGF5b3V0JHksIGNsYXlvdXQkY3QsIG1lYW4pKQogIGxhYmRmJGxhYiA9IHJvd25hbWVzKGxhYmRmKQogIGNvbmR4ID0gbGFiZGYkeTxtaW4odGVzdEB0cmVlJHNlZ21lbnQucHNldWRvdGltZS5saW1pdHMkZW5kKQogIGxhYmRmWyFjb25keCwieCJdID0gcm91bmQobGFiZGZbIWNvbmR4LCJ4Il0sIDApCiAgbGFiZGZbY29uZHgsIngiXSA9IDAuNQogIAogIHBsdDMgPSBwbHQyK3NjYWxlX2NvbG91cl9tYW51YWwodmFsdWVzID0gY29sc19jYywgbGltaXRzID0gZm9yY2UpKwogICAgZ2dyZXBlbDo6Z2VvbV9sYWJlbF9yZXBlbChkYXRhID0gbGFiZGYsIG1hcHBpbmcgPSBhZXMoeCA9IHgsIHkgPSB5LCBsYWJlbCA9IGxhYiksIHNpemUgPSAyLjIsIAogICAgICAgICAgICAgICBsYWJlbC5wYWRkaW5nID0gdW5pdCgwLjEsICJsaW5lcyIpLCBib3gucGFkZGluZyA9IHVuaXQoMC4xLCAibGluZXMiKSkrCiAgICBsYWJzKHRpdGxlID0gbikrCiAgICB0aGVtZShsZWdlbmQudGV4dCA9IGVsZW1lbnRfdGV4dChzaXplID0gNiksCiAgICAgICAgICBsZWdlbmQua2V5LnNpemUgPSB1bml0KDAuMiwgImxpbmVzIiksCiAgICAgICAgICBwbG90LnRpdGxlID0gZWxlbWVudF90ZXh0KHNpemUgPSA4KSkKICBwcmludChwbHQzKQogIHRyZWVfcGx0X2xbW25dXSA9IHBsdDMKfQpgYGAKClNhdmUgdHJlZXMKCmBgYHtyfQpmb3IobiBpbiBuYW1lcyh0cmVlX3BsdF9sKSl7CiAgcGRmKHBhc3RlMCgicmVzdWx0cy9wc2V1ZG90aW1lL0Rpdl8iLCBuLCAiX3RyZWUucGRmIiksIGhlaWdodCA9IDMuNSwgd2lkdGggPSA0KQogIHByaW50KHRyZWVfcGx0X2xbW25dXSkKICBkZXYub2ZmKCkKfQpgYGAKCkxvYWQgUk5BIHZlbG9jaXR5IHBzZXVkb3RpbWUKCmBgYHtyfQpnbHV0X2RhdF9kZiA9IHJlYWQuY3N2KCJyZXN1bHRzL1JOQXZlbG9jaXR5X0Rpdi9oZWF0bWFwX2dlbi9kaXZfZGF0X2RmLmNzdiIsIAogICAgICAgICAgICAgICAgICAgICAgIGhlYWRlciA9IFQsIHJvdy5uYW1lcyA9IDEpCm5ld2NlbGxuYW1lcyA9IHJvd25hbWVzKGdsdXRfZGF0X2RmKQpuZXdjZWxsbmFtZXMgPSBnc3ViKCJfMl93cGkiLCAiLTFfMSIsIG5ld2NlbGxuYW1lcykKbmV3Y2VsbG5hbWVzID0gZ3N1YigiXzRfd3BpIiwgIi0xXzIiLCBuZXdjZWxsbmFtZXMpCm5ld2NlbGxuYW1lcyA9IGdzdWIoIl82X3dwaSIsICItMV8zIiwgbmV3Y2VsbG5hbWVzKQpuZXdjZWxsbmFtZXMgPSBnc3ViKCJfOF93cGkiLCAiLTFfNCIsIG5ld2NlbGxuYW1lcykKbmV3Y2VsbG5hbWVzID0gZ3N1YigiXzEyX3dwaSIsICItMV81IiwgbmV3Y2VsbG5hbWVzKQpuZXdjZWxsbmFtZXMgPSBnc3ViKCJfMV93cGlfcG9zIiwgIi0xIiwgbmV3Y2VsbG5hbWVzKQpyb3duYW1lcyhnbHV0X2RhdF9kZikgPSBuZXdjZWxsbmFtZXMKYGBgCgpQbG90IHBzZXVkb3RpbWUgY29tcGFyaXNvbnMKCmBgYHtyfQpwbHRfY29tcF9sID0gbGlzdCgpCmZvcihuIGluIG5hbWVzKHN1YmZpeF9sKSl7CiAgcGxvdF9kZiA9IG1lcmdlKGdsdXRfZGF0X2RmLCBzdWJmaXhfNF9sW1tuXV1AcHNldWRvdGltZSwgYnkgPSAwKQogIGNwZSA9IHJvdW5kKGNvcihwbG90X2RmJG5ld3B0LCBwbG90X2RmJHBzZXVkb3RpbWVfNCwgbWV0aG9kID0gInBlIiksIDIpCiAgCiAgcGx0X2NvbXBfbFtbbl1dID0gZ2dwbG90KHBsb3RfZGYsIGFlcyh4ID0gbmV3cHQsIHkgPSBwc2V1ZG90aW1lXzQsIAogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgY29sb3VyID0gcHJlZF9jdGFsbCkpKwogICAgZ2VvbV9wb2ludCgpKwogICAgc2NhbGVfY29sb3VyX21hbnVhbCh2YWx1ZXMgPSBjb2xzX2NjLCBsaW1pdHMgPSBmb3JjZSkrCiAgICAjc2NhbGVfc2l6ZV9tYW51YWwodmFsdWVzID0gYygwLjUsIDAuOTUsIDEuMzUpKSsKICAgIGxhYnModGl0bGUgPSBuLCBzdWJ0aXRsZSA9IHBhc3RlMCgiUENDID0gIiwgY3BlKSwgCiAgICAgICAgIHggPSAiUHNldWRvdGltZSAoc2NWZWxvKSIsIHkgPSAiUHNldWRvdGltZSAoVVJEKSIpKwogICAgdGhlbWVfY2xhc3NpYygpKwogICAgdGhlbWUoYXNwZWN0LnJhdGlvID0gMSwKICAgICAgICAgIGF4aXMudGV4dCA9IGVsZW1lbnRfdGV4dChjb2xvdXIgPSAiYmxhY2siLCBzaXplID0gNiksCiAgICAgICAgICBheGlzLnRpdGxlID0gZWxlbWVudF90ZXh0KHNpemUgPSA3KSwKICAgICAgICAgIHBsb3QudGl0bGUgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDgpLAogICAgICAgICAgcGxvdC5zdWJ0aXRsZSA9IGVsZW1lbnRfdGV4dChzaXplID0gNyksCiAgICAgICAgICBsZWdlbmQudGV4dCA9IGVsZW1lbnRfdGV4dChzaXplID0gNiksCiAgICAgICAgICBsZWdlbmQudGl0bGUgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDcpLAogICAgICAgICAgbGVnZW5kLmtleS5zaXplID0gdW5pdCgwLjM1LCAiY20iKSkKfQpgYGAKClNhdmUgY29tcGFyaXNvbnMKCmBgYHtyfQpmb3IobiBpbiBuYW1lcyhwbHRfY29tcF9sKSl7CiAgcGRmKHBhc3RlMCgicmVzdWx0cy9wc2V1ZG90aW1lL0Rpdl8iLCBuLCAiX2NvbXAucGRmIiksIGhlaWdodCA9IDQsIHdpZHRoID0gNC41KQogIHByaW50KHBsdF9jb21wX2xbW25dXSkKICBkZXYub2ZmKCkKfQpgYGAKCk1ha2UgcHJvcG9ydGlvbiBwbG90cwoKYGBge3J9CnNtb19wcm9wX2xpc3QgPSBsaXN0KCkKZm9yKG4gaW4gbmFtZXMoc3ViZml4XzRfbCkpewogICMgc3Vic2V0IGRhdGEKICBzdWJtZXRhID0gZGF0YS5mcmFtZSgicHNldWRvdGltZSIgPSBzdWJmaXhfNF9sW1tuXV1AcHNldWRvdGltZSRwc2V1ZG90aW1lXzQsCiAgICAgICAgICAgICAgICAgICAgICAgImNlbGxjbHVzdGVycyIgPSBzdWJmaXhfNF9sW1tuXV1AbWV0YSRjZWxsY2x1c3RlcnMpCiAgCiAgbWVkX3B0X2NjID0gc29ydCh0YXBwbHkoc3VibWV0YSRwc2V1ZG90aW1lLCBzdWJtZXRhJGNlbGxjbHVzdGVycywgbWVkaWFuKSkKICAKICBsdF9iaW5zID0gY3V0KHN1Ym1ldGEkcHNldWRvdGltZSwgMTAwKSAjIDEwMCBlcXVhbGx5LXNpemVkIGJpbnMKICBwbG90X2RmID0gZGF0YS5mcmFtZShiaW5zID0gbHRfYmlucywgCiAgICAgICAgICAgICAgICAgICAgICAgY3N0ID0gYXMuY2hhcmFjdGVyKHN1Ym1ldGEkY2VsbGNsdXN0ZXJzKSkKICB0YWJfZGYgPSB0YWJsZShwbG90X2RmJGJpbnMsIHBsb3RfZGYkY3N0KQogIAogICMgcmVtb3ZlIGNlbGwgdHlwZXMgdGhhdCBhcmUgdG9vIHJhcmUgKDw1JSkKICB0YWJfZGYgPSByZXNoYXBlMjo6bWVsdCh0YWJfZGYvcm93U3Vtcyh0YWJfZGYpKQogIHVzZWNsID0gdGFwcGx5KHRhYl9kZiR2YWx1ZSwgdGFiX2RmJFZhcjIsIGZ1bmN0aW9uKHgpIGFueSh4PjAuMDUpKQogIHBsb3RfZGYgPSBwbG90X2RmW3Bsb3RfZGYkY3N0ICVpbiUgbmFtZXModXNlY2wpW3VzZWNsXSxdCiAgdGFiX2RmID0gdGFibGUocGxvdF9kZiRiaW5zLHBsb3RfZGYkY3N0KQogIAogICMgbm9ybWFsaXNlIGJ5IGNlbGwgdHlwZSBhYnVuZGFuY2UKICBwcm9wX3cgPSBwcm9wLnRhYmxlKHRhYmxlKHBsb3RfZGYkY3N0KSkKICB0YWJfZGYgPSB0KGFwcGx5KHRhYl9kZiwgMSwgZnVuY3Rpb24oeCkgeC9wcm9wX3dbY29sbmFtZXModGFiX2RmKV0pKQogIAogICMgcmVzaGFwZQogIHRhYl9kZiA9IHJlc2hhcGUyOjptZWx0KHRhYl9kZi9yb3dTdW1zKHRhYl9kZikpCiAgdGFiX2RmJFZhcjIgPSBhcy5jaGFyYWN0ZXIodGFiX2RmJFZhcjIpCiAgCiAgIyBwcmV2ZW50IGRpc2NvbnRpbnVpdHkgYnkgY29weWluZyB0aGUgcHJldmlvdXMgY29sdW1uIChsaWtlbHkgbm90IGhhcHBlbmluZykKICB0YWJfZGYgPSB0YWJfZGZbb3JkZXIodGFiX2RmJFZhcjEsIGRlY3JlYXNpbmcgPSBGKSxdCiAgZm9yKGkgaW4gdW5pcXVlKHRhYl9kZiRWYXIxKSl7CiAgICBpZihhbnkoaXMubmFuKHRhYl9kZiR2YWx1ZVt0YWJfZGYkVmFyMT09aV0pKSl7CiAgICAgIHRhYl9kZiR2YWx1ZVt0YWJfZGYkVmFyMT09aV0gPSBwcmV2CiAgICB9CiAgICBwcmV2ID0gdGFiX2RmJHZhbHVlW3RhYl9kZiRWYXIxPT1pXQogIH0KICAKICAjIHNtb290aGVuIHRoZSBwcm9wb3J0aW9ucyAoYW5kIGZvcmNlIGNvbnN0cmFpbiB0byAwLTEpCiAgdGFiX2RmMiA9IHRhYl9kZgogIHRhYl9kZjIkdmFsdWUyID0gdGFiX2RmMiR2YWx1ZQogIGZvcihpIGluIHVuaXF1ZSh0YWJfZGYyJFZhcjIpKXsKICAgIGZmZiA9IGxvZXNzKHZhbHVlfmFzLm51bWVyaWMoVmFyMSksIGRhdGEgPSB0YWJfZGYyW3RhYl9kZjIkVmFyMj09aSxdLCAKICAgICAgICAgICAgICAgIHNwYW4gPSAwLjUpCiAgICBwcmVkID0gcHJlZGljdChmZmYpCiAgICBwcmVkW3ByZWQ+MV0gPSAxCiAgICBwcmVkW3ByZWQ8MF0gPSAwCiAgICB0YWJfZGYyJHZhbHVlMlt0YWJfZGYyJFZhcjI9PWldID0gcHJlZAogIH0KICAKICAjIGZvcmNlIGNvbnN0cmFpbiBlYWNoIGludGVydmFsIHRvIDAtMSBieSBkb2luZyBwcm9wb3J0aW9uCiAgZm9yKGkgaW4gdW5pcXVlKHRhYl9kZjIkVmFyMSkpewogICAgdGFiX2RmMiR2YWx1ZTJbdGFiX2RmMiRWYXIxPT1pXSA9IHRhYl9kZjIkdmFsdWUyW3RhYl9kZjIkVmFyMT09aV0vc3VtKHRhYl9kZjIkdmFsdWUyW3RhYl9kZjIkVmFyMT09aV0pCiAgfQogIAogIHRhYl9kZjIkbWFqb3IgPSB1bmxpc3QobGFwcGx5KHN0cnNwbGl0KHRhYl9kZjIkVmFyMiwgIl8iKSwgZnVuY3Rpb24oeCkgeFsxXSkpCiAgcmVzID0gbGlzdCgpCiAgZm9yKG5ubiBpbiB1bmlxdWUodGFiX2RmMiRWYXIxKSl7CiAgICBzcz10YXBwbHkodGFiX2RmMlt0YWJfZGYyJFZhcjE9PW5ubiwidmFsdWUyIl0sIHRhYl9kZjJbdGFiX2RmMiRWYXIxPT1ubm4sIm1ham9yIl0sIHN1bSkKICAgIHJlc1tbbm5uXV0gPSB3aGljaC5tYXgoc3MpCiAgfQoKICBzbW9fcHJvcF9saXN0W1tuXV0gPSBnZ3Bsb3QodGFiX2RmMiwgYWVzKHggPSBWYXIxLCB5ID0gdmFsdWUyLCBncm91cCA9IFZhcjIsIGZpbGwgPSBWYXIyKSkrCiAgICBnZW9tX2FyZWEoKSsKICAgIHNjYWxlX3lfY29udGludW91cyhleHBhbmQgPSBjKDAsMCkpKwogICAgc2NhbGVfZmlsbF9tYW51YWwodmFsdWVzID0gY29sc19jY1tuYW1lcyhjb2xzX2NjKSAlaW4lIHRhYl9kZjIkVmFyMl0pKwogICAgbGFicyh4ID0gIkJpbnMiLCB5ID0gIlByb3BvcnRpb24iLCBmaWxsID0gIkNlbGwgdHlwZSIpKwogICAgdGhlbWVfY2xhc3NpYygpKwogICAgdGhlbWUoYXhpcy50ZXh0LnggPSBlbGVtZW50X2JsYW5rKCksCiAgICAgICAgICBheGlzLnRleHQueSA9IGVsZW1lbnRfdGV4dChzaXplID0gNi41LCBjb2xvdXIgPSAiYmxhY2siKSwKICAgICAgICAgIGF4aXMudGlja3MueCA9IGVsZW1lbnRfYmxhbmsoKSwKICAgICAgICAgIGF4aXMubGluZSA9IGVsZW1lbnRfYmxhbmsoKSwKICAgICAgICAgIGF4aXMudGl0bGUgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDcpLAogICAgICAgICAgbGVnZW5kLnRleHQgPSBlbGVtZW50X3RleHQoc2l6ZSA9IDYpLAogICAgICAgICAgbGVnZW5kLnRpdGxlID0gZWxlbWVudF90ZXh0KHNpemUgPSA3KSwKICAgICAgICAgIGxlZ2VuZC5rZXkuc2l6ZSA9IHVuaXQoMC40LCAiY20iKSkKfQpgYGAKClNhdmUgcHJvcG9ydGlvbnMKCmBgYHtyfQpmb3IobiBpbiBuYW1lcyhzbW9fcHJvcF9saXN0KSl7CiAgcGRmKHBhc3RlMCgicmVzdWx0cy9wc2V1ZG90aW1lL0Rpdl8iLCBuLCAiX3Byb3AucGRmIiksIGhlaWdodCA9IDMsIHdpZHRoID0gNS4yKQogIHByaW50KHNtb19wcm9wX2xpc3RbW25dXSkKICBkZXYub2ZmKCkKfQpgYGAK